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DEFINITION OF SYMBOLS 


Symbol 

A 

ai 2 , a 2 2 
bj 2 , b 2 2 

c l> c 2 > c 3 

E 

E 

e 

F 

f 

8 

H 

^1* ^2» J 3» ^4 
K 1# K 2 
L 


Definition 

Constant measuring force of planar field 

Larger roots of P(£) and P(t?), respectively 

Smaller roots of P(£ ) and P(r )) , respectively 

Integration constants used in solution of Hamilton- 
Jacobi equation 

Separation constant in Hamilton- Jacobi equation 

When used with arguments, denotes an elliptic 
integral of the second kind 

Denotes eccentricity of an unperturbed Keplerian 
orbit 

Constant measuring force of planar field; F = mA 
Any continuous function of £ 

Any continuous function of r\ 

Hamiltonian 

Integrals defined by equations (69) , (71) , (74) , and 
(76) 

Parameters in elliptic integrals associated, 
respectively, with £ and 17 . 

Vectorial angular momentum in an unperturbed 
Keplerian orbit 
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DEFINITION OF SYMBOLS (Continued) 


Symbol 

L 


1 


m 


P(S) 

P(t?) 


q i 


r 

r 

S 

S l> S 2 
t 

u 

u 


Definition 

Magnitude of angular momentum in an unperturbed 
Keplerian orbit 

Separation constant in the Hamilton-Jacobi equation 
Mass of satellite 


Denotes polynomial in specifically, as used, a 
quadratic polynomial 

Denotes polynomial in rj: specifically, as used, a 
quadratic polynomial 

Denotes any of the momenta p , p , p defined by 


equations ( 12) , ( 13) , and ( 14) 




Denotes any of the momenta conjugate to p , P^, and 

P 0 


Radius vector from primary to satellite 
Magnitude of radius vector from primary to satellite 
Jacobi substitution function 

Components of S which are functions only of £ and rj, 
respectively 

Time, a parameter 

Potential function 

Dummy variable used in evaluation of I I( . . ., 14 ; 
different substitutions are used for each integral 


Symbol 

W 


DEFINITION OF SYMBOLS (Continued) 

Definition 

Laplace integration constant for unperturbed two-body 
problem 


W Magnitude of integration constant for unperturbed 

two-body problem 

W Function related to S by equation ( 24) 

x,y,z Cartesian coordinates describing motion of satellite 

x,y Dummy variables 


/ 3 Separation constant in Hamilton- Jacobi equation 

rj One of three parabolic coordinates, rj = r - z 

0 Polar angle in spherical coordinates 

H Gravitational parameter of attracting primary 

£ One of three parabolic coordinates £ = r + z 

p Projection of radius vector in (x, y) plane 

< p Angular coordinate from arbitrary zero point in 

(x, y) plane 

ip Angle between Laplace vector integration constant and 

the radius vector 
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DEFINITION OF SYMBOLS (Concluded) 


Symbol 

Subscripts 

i 

0 

x,y,z 

Other Notation 


x 


Definition 

Dummy subscript or summation variable 
Indicates initial point 

Indicates Cartesian component of a vector 

Indicates derivative with respect to time 
Indicates a scalar product 
Indicates a vector product 
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MOTION OF AN ARTIFICIAL SATELLITE UNDER 
COMBINED INFLUENCE OF PLANAR AND 
KEPLER IAN FORCE FIELDS 


INTRODUCTION 


The development of the technology of artificial earth satellites has 
posed mathematical problems which heretofore were academic. One such 
problem is the radiation pressure perturbations of satellites with large area- 
to-mass ratios (i. e. , Echo) . The standard method of treatment of such a 
problem is a perturbational analysis of a purely Keplerian orbit which is dis- 
turbed by a radiation pressure field. 

It will be shown in this report that a generalized form of the two-body 
problem which more closely approximates the radiation pressure problem can 
be solved in closed form. Specifically, the problem of the motion of a satellite 
under the combined influence of Keplerian and planar force fields will be 
solved in terms of the elliptic functions. 

The assumption of a planar force field is still an approximation to the 
true radiation pressure problem ’’ because the divergence of the field is 
neglected. It is obvious that neglecting the divergence of the radiation pressure 
field is far less restrictive than disregarding the entire field. Perturbational 
treatment with the present solution as a model should converge far more 
readily than treatments which assume an initially Keplerian model. 


'r* 

The treatment here has equal application to the case of the classical Stark 
effect; that is, the radiation from hydrogen gas in an external planar electric 
field is (classically) analogous to the radiation pressure problem. See 
E. Schrodinger, Ann. Physik, 80, 457, 1926, and P. S. Epstein, Phys. Rev., 
28, 695, 1926. 


The solution may be outlined as follows. The equation of motion in 
Cartesian coordinates is used to formulate the Hamilton- Jacobi equation of the 
problem. This equation is then rewritten in terms of parabolic coordinates 
and separated by means of the usual product of independent functions assumption.* 
The solution is effected by quadratures of the separated equations in the three- 
dimensional case. These quadratures are explicitly integrated for the case of 
two-dimensional motion, which is initially circular. A graphical presentation 
of the results is given. 


FORMULATION OF PROBLEM 


The equations of motion of a satellite under the combined influence of 
an inverse square field and a uniform planar field are 


.. _ MX 

X " '(x 2 + y 2 + z 2 ) 3/2 

MV , 

y - - (x 2 + y2 + z 2 )3 A > 

and 

z - - (x 2 +y 2 +z 2 )3 /2 +A . 

where x,y,z are Cartesian coordinates, m the gravitational parameter of the 
attracting center, and A the planar force-to-mass ratio. The force may be 
derived from a potential U defined by 

U= - U(x,y, z, A) = -m ( (x 2 + ^" 2 ^ 1/2 + Az ) . 


The problem may now be stated as follows. Determine functional 
relationships which relate x,y, and z to the time, initial conditions, and param- 
eters M and A. This procedure is quite difficult, if not impossible, in Cartesian 
coordinates. For this reason, we change from Cartesian coordinates to 
parabolic coordinates. 


The separation of the Hamilton- Jacobi equation in these coordinates is 
suggested in Reference 1. It should be noted that the problem of a planar field 
superimposed on an inverse square field is a special case of Euler s problem 
of two fixed centers. See Euler, Mem. de Berlin, 1760, p. 228. 
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TRANSFORMATION OF COORDINATES AND THE 
HAMILTON-JACOBI EQUATION 


This potential may be rewritten in cylindrical coordinates as 


U = U( p, 0 , z) = Z— - mAz 


(i) 


where 


and 


x = p cos 
y = p sin cp 
z = z 


(2) 


The Lagrangian yielding the basic set of differential equations may be 
written as 


m , »o p*o .o m/i 

L = — (p 2 + p 2 <p 2 + z 2 ) + — 

(p 2 + z 2 ) 2 

= (p 2 + P 2 </> 2 + z 2 ) - U(p,<p, zj * 


+ mAz 


( 3 ) 


The Lagrangian is now transformed to parabolic coordinates £, p,0. 
This will ultimately allow a separation of the Hamilton-Jacobi equation 

p = -/Jr) 

„ = 4- n 


and 


<P = <P 

The radius vector in spherical coordinates, r, is given by 


( 4 ) 

( 5 ) 

( 6 ) 


■ J p 2 +2 2 


( 7 ) 


* The inclusion of m as a multiplicative factor on A is, of course, unnecessary. 
It simply turns out to be algebraically convenient. It can be removed in sub- 
sequent equations by setting A = F/m if desired. 
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Thus, 


and 


2 \ nprf 


( 8 ) 


z = 


JLzA 


(9) 


so that 


mU+h) /if + + m <j> 


8 






- U(^,rj,0) . 


( 10 ) 


Defining the momenta conjugate to rj, 0 as 
3 L 

P i 9q 4 


(li) 


(where q. is of any of rj, or (p and p. is any of p , p , or p Iwe find that 


r 




and 


D _ m(^) 


m(^+rj) • 

p = — ^ rj , 
7] 4rj 


P^ = m^r;0 


so that 


4 £ P 


£ = 


L 


m( |+rj) 


(12) 

(13) 

(14) 

(15) 


and 


Now 


V = 


<P = 


411 p „ 

m( £+ 77 ) 


jL 


m £ tj 


I’W >7 p ^ 2 > 


L(|. V, <t>, P t , P„, Pj - 


(16) 


(17) 


2m | 7 } 


- U(£, 7 ?, 0) . 


(18) 
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The Hamiltonian, H, is defined by 


2 


H = 


T q. p. 
. i i 
l 


L(q,p) = 


P i 


m (I+17) 


Up ^ 2 


+ V P^ 2 ) 


2m |r? 


+ U(|,tj,0) . (19) 

We now define a function S such that 



p r 

as 

(20) 


p , = 

as 

dv 

(21) 

and 

p <r 

as 

a^ 

(22) 


The solutions to the original problem now must be solutions to the 
Hamilton- Jacobi equation 


rnU+T]) 


Kf-)-(fr) 2 ] 


2 m £ 17 \d<p / 


(■fry + 


9S 

at 


(23) 


SEPARATION OF HAMILTON- JACOB I EQUATION: 
SOLUTION OF PROBLEM 


We shall attempt to find a solution to equation (23) in the form 


S(£,tj,<M) = W( 4,rj,0) - Et 


(24) 


where E is a constant. Then 



2 


m (£+r?) 



1 

2 m £77 





= E . 


(25) 


It is not necessary, as yet, to specify the exact form of U(|,tj, 0 ) . 
Assume, more generally, that 


U(£,rj,0) 


f(£ ) + g(T? ) 
£ + TJ 


(26) 


where f(£) and g(rj) are simply well-behaved functions of their arguments which 
allow the integrals to exist. Inserting ( 26) into ( 25) and setting 


gives 


W = 1 4> + Sid) + S 2 (tj) 


2 

t /9SA 2 , „ /9S, V 

, l 2 , f(i ) + R(ri) 

m( £+ 77 ) 

* \ 9£ / ^ V 9t } / 

2 m £ 77 £ + 77 


(27) 


= E , (28) 


Multiplying (28) by md+ 77 ) and rearranging gives 


25 (if ) 2 + mf(5) - mE « + if = ' z ” 


+ m E rj - 


2r) 


(29) 


By the usual argument employed at this state in the solution of a partial 
differential equation, we note that the left side of (29) depends only upon | , 
while the right side depends only upon rj. Since £ and 77 are independent coordi- 
nates, we can move £ through an arbitrary sequence of values while 77 is constant, 
then reverse the process. But the equality given by (29) must be retained so 
that both sides must be equal to a constant. 
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We, thus, have 


and 


2£ 


2r) 


2£ 


'dS2 X “ 


dr] 


+ m g(rj) - m Er] + 


27 ] 


= - (3 


(30) 


(31) 


Equations (30) and (31) integrate to yield 


and 


S i = / 
lo 

v 

s, = / 

Vo 




2rj 


m f (£) 

l 2 


2£ 

4| 2 

0 

mg(rj) 

l 2 

drj 


d« 


2r] 


4rf 


(32) 


(33) 


Inserting (32) and (33) into (27) and the resultant equation into (24) gives 

J 1 J HlE + _ mf({) _ _j! 


S(4,r],^,t) = 10 - Et + 


!( 


2£ 


2^ 


4^ 




V 

+ I 

Vo 


m E / 3 m g(r?) I 2 


dr/ 


(34) 


2r/ 2 tj 40 

We now specialize f(£) and g(r/ ) by writing equation ( 1) in the form 

(-n* “0)_( m ,^) 


-m p 

2.U r,2 


\/p 2 + z 


mAz = 


-2mp 

4+r] 


m A ( g -T] ) 
2 


I + 1? 


( 35 ) 
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Comparison of equations (26) and (35) shows that 


f(|) - - mu - — 

(36) 

and 


, . mAr) 2 

g(rj) mp + 2 

(37) 


so that (34) becomes 


_ r’ I m E (13 + mV) . m 2 A£ 

s(^nAt) =i*-Et+ J J — + ^ E + £ - 


2£ 


4r 


d£ 


V 

+ I „ 

n o 


m E ( -P + m 2 ju) 
2 277 


m 


2 Ar) 


4rj 2 


dr? 


(38) 


The arbitrary constants in this equation are 1, E, and (3. Following the 
normal procedures of Hamilton- Jacobi mechanics, we now equate the partial 
derivatives of S with respect to 1 or E or (3 to additional constants yielding the 
three equations 


, , m r 

'1 = ~ t+ ^ J 


AA L 


2 £ 0 '\/+m 2 A£ 3 + 2mE£ 2 + 2(/?+m 2 |tt)£ - l 2 


7) dr) 


m j 

2 Tj 0 \l -m 2 A?) 3 + 2mEr) 2 + 2(m 2 p-/?) rj - l 2 


(39) 


- / 


dl 


£ 0 4+m 2 A^ + 2mE£ 2 + 2(j3+m 2 p) £ - l 2 


-/ 


dr) 




\/- m 2 A?) 3 + 2mEt] 3 + 2(m 2 p-/3)r? - l 2 


(40) 
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and 







t.V + m 2 A | 8 + 2mE | 2 + 2(/3 + m 2 p) £ - l 2 


dr] 


i 

f - ... 

T}t r] 'sI-m 2 Ari z + 2mErj 2 + 2(m 2 77 -/ 3)77 - l 2 ( ’ (41) 


Choosing £ = rj =.7]$, we find that Cj is the initial time t 0 . Similarly, c 2 
vanishes because we initialize motion at a possible orbital point (£ 0 ,t 7 0 ). 

Finally, c 3 is the initial value of <p. Equation (39) can then be regarded as 
the relation maintained between | and 77 throughout the motion. If this equation 
is solved to yield 77 = 77 (ij), we have the projection of the motion in the |-tj plane. 
Then equation (41) is the full three-dimensional description of the spatial motion 
which relates 4> to £ and 77 . Elimination of 17 ( or £) via equation ( 40) would yield 
0(£) (or 0 ( 77 )) . Finally, equation (39) predicts values of the parameter t for 
known values of £ and r\. We can thus obtain values for any three members of 
the set {|, 77 , 0 ,t} once any one of the four is specified. 


It should be noted that, in a computational sense, we can obtain answers 
much more easily if a value of 4 or 77 is originally specified rather than a value 
of <p or t. This observation is not apparent from the material developed thus 
far but will be indicated by a special case to be treated in detail below. It 
arises because of the difficulty of inverting the transcendental expressions which 
result from the integration of the quadratures appearing in equations ( 39 ) and 
(41). 


We would expect equation (39) to be complicated because time is usually 
difficult to treat in celestial mechanics problems. Equation (41) relates all 
three spatial parameters and, for this reason, may be expected to be more 
cumbersome than a relationship between two spatial parameters, equation (40) . 

The next problem of interest is the physical identification of the constants, 
1, E, and /3. 
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IDENTIFICATION OF CONSTANTS IN TERMS OF 
INITIAL CONDITIONS 


The constant E, first introduced in equation (24) , can be most easily 
identified in Cartesian coordinates. Since we have 


X = \l £, TJ cos <fi 
y - n/T V sin <p 

and 



equation (25) becomes 


(42) 

(43) 

(44) 


E = 


m(£ +r]) 

m(£+n) ( JL 


7 P^J + 

JL 

2m^r] 

jiL) + 

m 

rj ' 

2 


8 ' $ 


( ^ + y2 + z 2 ) - ^ 


P 2 + U,T? .0) 


+ UU,TJ,0) 


- mAz 


(45) 


^ N/ x 2 +y 2 +z 2 


Notice that the first two terms following that last equality in equation (45) represent 
the energy in an unperturbed two-body problem. 


For the identification of 1 we turn to equations (22) , (24) , and (27) to 
find that 



Applying equations (42) and (43) gives 


1 = m(x 2 + y 2 ) 4> 


m(xy - yx) 


( 47 ) 


so that 1 is simply a component of angular momentum. 

The constant requires more work for its identification. Note that from 


and 


so that 


(20) and (21) 

/ as \ 2 

mE 

(/S+m 2 p) 

+ m2A 4 

l 2 

(48) 

v as / 

2 

24 

4 

44 2 

/ as\ 2 

mE 

( -/3+m 2 p) 

m 2 A r] 

l 2 

(49) 

V arj / 

2 

277 

4 

4t) 2 

2 _ JL 

/4 + 7\ 

m 2 p / 

4 -7 \ 

ItlHl 


- ? 4 “ ' 2 

V 47 / 

2 * 

‘ 47 ' ' 

4 ( £ + T ?> 



I 2 - n 2 
4 V 


(50) 


Solving for /3 gives 


P ~ + 


247 
£ + 7 


( p /- p v ) 


p, ! - p-M + (f^) 


1 2 U -rj) ( 
47 j 


m 2 A£r) 


(51) 


This appears much less elegant than the expressions derived for E and 1 
so that a further comment is in order. In Appendix A it is shown that /3 is a 
modification of one component of Laplace's vectorial integration constant. 
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We have now arrived at such a point in the development that no further 
progress can be made without specializing the general results derived thus 
far. This is because the integrated forms of equations (39) , (40) , and (41) 
cannot be obtained until such time as certain restrictions on the magnitude of A 
and the initial conditions are spelled out. For any given set of such conditions, 
the integration can proceed but the conditions must be given. 

To illustrate the procedure, we turn our attention to the special case 
of two-dimensional motion with initially circular orbit conditions. 


INTEGRATION OF EQUATIONS OF 
TWO-DIMENSIONAL MOTION 


The choice of initial coordinates was made in such a way that the planar 
force acts in the direction of increasing z. We can then expect that if motion 
is to occur in a plane, this plane must include the z axis. It is algebraically 
simplest to define the plane of two-dimensional motion to be the x-z plane. 

This, however, is simply a convenience. If such a choice is made, then the 
coordinate <p is identically zero. If <f> = 0, then certainly 0 = 0, and equation 

Using this condition, equations (39) and (40) become 


d 4 r . _ 

N /m 2 A| 3 + 2mE£ + 2(/3 + m 2 | u) 


(46) shows that 1 = 0. 




f 




and 



/ 

Vo 


\frj dri __ 

\J -w 2 Ar) 2 + 2mEr; + 2{m 2 jj.-(3) 


(52) 


/ 


jiL 


| 0 m 2 A| 2 + 2mE£ + 2(/3+m 2 //)] 


dr? 


= / 

tj 0 Jr? [-m 2 Ar7 2 + 2mEr? + 2(m 2 ?u-?S)] 


(53) 
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To simplify the algebra, we now make use of the fact that the choice of 
units for mass, time, and length are at our disposal. Choose these units in 
such a way that 


m = 1 

M = X 

r 0 = 1 (54) 

Having made this choice of units, we further choose initial conditions which 
ensure an (initially) circular orbit of radius r 0 . Then 


x 0 = 1 

Yo = 0 

zo = 0 

x 0 = 0 (55) 

Yo = o 

z 0 = ±1 

where a subscript of zero indicates an initial condition. Using these conditions 
in equations (45) and (51)*, we find 

J_ 

2 

A_ 

2 (56) 



First, convert initial conditions to ^ 0 , 77 0 coordinates. 
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Using equations ( 54) and ( 56) , the polynomials appearing in the 
denominator of equations ( 52) and ( 53) may be written as 


P(£) = m 2 A £ 2 + 2 mE£ + 2(/?+m 2 M) = A £ 2 



A_ 

2 


) 


= A ( a t - |) (bj - |) (57) 

and 

P(tj) = -m 2 Arj 2 + 2 mE ?7 + 2 (m 2 ^-/?) = - A 77 2 - 77 + 2 ^1 + 

= A(a 2 - ??) ( | b 2 | + 77 ) (58) 


where a l5 bj and a 2 , b 2 are the roots of P(£)and P(r/), respectively. 
Equations ( 52) and ( 53) may now be written as 


t - t n = 


and 


/ 




/ 


\frf di] 


^0 ^/Ta7^TMb7^T)” ^/( a 2 - 77 ) ( |b 2 | + 77 ) 


(59) 


/ 


.di 


4o ^£( a l ~ £) - £) 


" / 


d7? 


(60) 


770 ^7 (a 2 - 77 ) ( | b 2 | + 77 ) 


We now examine the possible singularities of the integrands of the last 
two equations. The integrals will exist only under certain restrictions, and a 
detailed examination is in order. The equations which define cylindrical 
coordinates in terms of £ and 77 , equations (3) and (4) , give 


and 


£ = z ± *Jx 2 + y 2 + z 2 


77 = -z ± ^x 2 + y 2 + z 2 . 
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Either 4 — 0 and 77 ^ 0 or | s: o and i) ^ o as may be seen from equation 
(3). Choose, for convenience, 


£ = z + n/x 2 + y 2 + z 2 

and 

7 ] =-z+ \/x 2 + y 2 + z 2 . (61) 

Now, £= : 0ifx = y= : z = 0orifx = y= 0 and z < 0. Similarly, 77 will vanish for 
x=y = z= :: 0 orx=y -0 and z > 0 . This gives £ ^ 0 and 17 S: 0 which establishes 
a lower bound on £ and 77. 

Label the roots of P( |) in equation (57) in such a way that a t ^ bj. With 
this convention we can now establish upper bounds on both £ and 17 . To do this, 
we notice that from equations (30) , (31) , (57) , and (58) we have 



Using the non-negativity property of £ and 77 , we may conclude that 
p(£) 2 ; 0 (63) 

and 

p(r]) 2 : 0 ( 64) 


so that 


I ^ 


or 


4 a bj 
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indicating a forbidden range of £ values. That is, we may expect that the 
region 


a i < £ < b i (65) 

is forbidden. 

Similarly, the non-negativity of P (r?) ensures that 
V — a 2 • 

We have thus established the bounds on £ and 77 as 

0 < | £ aj ( 66) 

and 

0 < r) — a 2 • (67) 

It is now apparent that equations (59) and (60) may possess singularities 
in the integrand. This means that these integrals might necessarily be regarded 
as improper (Riemann) integrals, however, and the existence of these integrals 
may be established by standard techniques of elementary calculus. 

Next, we shall consider the problem of bounds on A. Let £* and 77 * be 
roots of P(|) and P(tj). Then 

PU*) = A** 2 - + 2^1 - y) 

and 

P( 77 *) - - At t * 2 - 77 * + 2 (l + . 

Then 4 * = £* (A), tj* = tj*(A). Specifically, 
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and 


rj* 


1 ± */ 1 - 8A (l - y) 
2A 

1 ± Jl + 8A (l + 

-2A 


It is not necessary to require that £* and 17* be real. The integrands, how- 
ever, change character if we allow imaginary values of the roots. For con- 
venience, therefore, we require that and 17* be real. Now 17* is real for 
A ^ 0. The. values of will be imaginary if 


1 


2 


< A < 1 + 


n/ T 

2 * 


Thus, we bound A within the limits 


0 < A < 1 


2 


For reasons of orientation, we include the following table. 


( 68 ) 


> 

a i 

bi 

a 2 

b 2 

-4 

3. 73 

3. 73 

1. 73 

-9. 19 

0. i 

2. 55 

7. 45 

1. 78 

-11. 83 

0. 01 

2. 03 

97. 97 

1. 97 

-101. 97 

0. 001 

2. 005 

998 

2 . 00 

-1002 

0. 0001 

2 . 000 

.9995 

2 . 00 

-10002 


( Note b 2 < 0 for all allowable values of A. ) 
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To reiterate the results of this section to this point, we now have 
established the inequalities 


and 


0 — | ^ a 4 
0 < 77 < a 2 
0 < A < 1 " 


» 


N (W 
2 


In equation ( 60) , we set 

h = / ^ (69) 

4o - £) (bj - O 


When we introduce a dummy variable u defined by 
| = a t sn 2 u 

with the elliptic parameter K given by 


K 2 = 


ai_ 

b i 


equation (69) becomes [3] 



Define 


h 



Vo 


drj 

n/ 17 ( a 2 - 17 ) ( |b 2 |+ v) 


Set 


(70) 


(71) 


77 = a 2 cn 2 u 
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with 


a 2 

a 2 + [bj| 


K 2 = 

Then, [ 3] 



Substituting equations (70) and (72) into equation (60) gives the 
equation of the orbit: 



Notice that we may not remove the elliptic functions in equation (60) by application 
of the operators sn or cn due to the multiplicative factors which occur. If 


a 2 


+ ft>2 


happens to be rational, this manipulation could be accomplished, if 


desired. 


Notice that neither transformation used in the derivation of equation (73) 
possesses a singularity. Equation (73) is a general solution valid for all 
£ and r). 
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Next, set 


I 3 



£0 


'll d£ 

J( aj - a (b ± - 0 


( 74 ) 


Setting 


£ = a t sn 2 u 

and 



we obtain [ 3] 



Finally, set 


1 

/ 

7o 


\frT dr; 


\l(a 2 - V) ( Ib 2 l + n) 


(75) 


(76) 


Setting 


77 = a 2 cn 2 u 


and 


K 2 


a 2 

a 2 + ! b 2 1 
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yields [ 3] 



Substituting equations (75) and (77) into equation (59) yields 



A comparison of equations (73) and (78) indicates that the problem of 
solving equation (73) for £ or rj is simpler than solving equation (78) for | or rj. 
It would thus be easier to write equations of the form 

t = t(£ ) 
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or 

t = t( 77 ) 
than to write 


£ = !<t) 


or 


V = v(t) . 

The difficulty with specification of £ or rj and, subsequently, solving 
for t is that we must choose a value of £ or 17 which exists for some value of t. 
In other words, only those values of | and rj which actually constitute a point 
on the orbit will yield a value of t. 


NUMERICAL RESULTS 


The preceding theory has dealt at some length with the solution of a 
seemingly trivial modification of the Kepler problem. This problem is defined 
by solution of the following differential equations. 


x = 


~JL 


(x 2 + y 2 + z 2 ) 


2v3/2 


(79) 


(x 


2 +y 2 + 


2 ^ 3/2 


Z*) 


( 80 ) 


z = 


( x 2 + y 2 + z 2 


)3/2 


+ A 


( 81 ) 


The results of a separation of variables of the last three equations is 
given by equations (39) , (40) , and (41) . For the special case of two-dimensional 
motion, which is initially a circular orbit, and A restricted to the range 


0 < A < 1- 


\HT 


the separated equations were integrated analytically. 
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The results obtained are fairly complicated, and it is worthwhile to obtain 
a graphical form of the solution. These graphical results are presented in this 
section. 

The following graphs were prepared by integration of equations (79) , 

(80) , and (81) . * For illustrative purposes a two-dimensional case was taken 
with initial conditions given by equation ( 55) . The values of A were chosen to 
give a representative set of orbits. Values of A were 

A = {0, 0.001, 0.01, 0.04, 0.05, 0.06, 0.08, 0.1} . 

For every value except the first two, two values of z (namely +1 and -1) were 
chosen. Furthermore, in every case, the results are plotted in both x-z and 
£-17 coordinate systems. 

The first set of graphs, Figures 1 and 2, is shown for the case of A = 0. 

In x-y coordinates, the result is a circle but in £-17 coordinates the result is a 
straight line. This illustrates the fact that results are strikingly different in 
the two-coordinate system. 

The second set of results. Figures 3 and 4, has an A value of 0. 001. 

The computer was stopped at an arbitrary point in the integration because of 
excessive running times (approximately 30 minutes on a SDS 930 computer 
utilizing a Cal Comp plotter) . The obvious characteristic of Figure 3 is that 
the orbit is shifting in the direction of -x, even though the planar force field is 
acting in the direction of +z. This condition illustrates the gyroscopic effect of 
orbital motion that is usually less apparent. In Figure 4, as in all other cases, 
the motion is limited by the roots of the polynomial, as was shown in equations 
(66) and (67) . 

A second important point in connection with the numerical integration is 
illustrated by Figures 5 and 6, for example. The motion, in general, tends 


'The reason for numerically integrating these equations instead of using equa- 
tions (73) and (78) was that existing numerical tables for elliptic functions are 
quite inadequate. A serious study of the behavior of the orbits under discussion 
would require series generation of the inverse elliptic functions. The above 
theory was checked by comparison of numerical integration of equations ( 52) and 
(53) with numerical integration of equations (79) and (81) (planar case). Equa- 
tions (73) and (78) were analytically differentiated to establish their validity. 
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to spiral the point mass toward the attracting center. The integration step size 
was reduced as the particle approached the center of attraction. The integration 
was terminated when the constants E or /? changed value across an integration 
step. (E was found to be more sensitive than /3 to such changes.) A careful 
examination of the behavior of the motion of the point mass in close proximity 
with the primary could best be examined via the analytical solution 

Without consideration of the remaining graphs (Figs. 7-28) in detail, 
we shall examine two additional points. The first of these is that in which the 
combination of forces is such that a backward loop may be attained by the 
particle. This effect is easily seen in Figure 11 by examining the point of 
greatest extension of motion in the x direction. 

Figures 27 and 28 illustrate a fairly interesting case. For an initially 
circular orbit, as always, and z = -1, an extremely close approach (0. 000005 
unit) occurs during the first orbit. The exact value of A for which a collision 
occurs during the first orbit could be determined by use of equation (73) if it 
is of interest. 
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FIGURE 5. Z 0 = 1, A = O.Oi CARTESIAN COORDINATES 



FIGUKE 6. Z 0 = 1, A = 0. 01 PARABOLIC COORDINATES 











“ FIGURE 13. Z 0 = 1, A = 0. 08 CARTESIAN COORDINATES 
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FIGURE 18. Z 0 = -1, A = 0. 01 PARABOLIC COORDINATES 
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FIGURE 21. Z 0 = -1, A = 0.05 CARTESIAN COORDINATES 
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FIGURE 22. Z 0 = -1, A = 0.05 PARABOLIC COORDINATES 
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FIGURE 24. Z 0 = -i, A = 0. 06 PARABOLIC COORDINATES 



08 CARTESIAN COORDINATES 






CONCLUSIONS 


It has been shown that orbital motion which involves combined Keplerian 
and planar force fields can be expressed analytically. The results are given 
in the form of quadratures for general three-dimensional motion. 

For the special case of two-dimensional motion and initially circular 
orbits these integrals are evaluated in terms of elliptic functions. Numerical 
results for these conditions are presented graphically. The varied results which 
are due to variations in the ratio of field strengths indicate the generalized 
form of the analytical results. 

The analytical results derived above were not used to generate the 
graphical results which are presented, due to the inadequacy of published 
numerical tables in the area of elliptic functions. However, the treatment 
of unusual behavior such as collisions and long term stability can be approached 
much more easily via the analytic solution. 

The above results are not necessarily directly applicable from the 
physical point of view. However, they do form a basis for further work in 
perturbation theory and for cases in which the combination of fields which 
is assumed in this report is a sufficiently accurate model. 

George C. Marshall Space Flight Center 

National Aeronautics and Space Administration 
Huntsville, Alabama, March 19, 1968 
125-17-05-00-62 


APPENDIX 


FURTHER CONSIDERATIONS OF SEPARATION CONSTANT (S 

Suppose that we consider the unperturbed two-body problem. Let the 
position be described by a radius vector "r . The equations of motion for this 
problem are 


r = - M 


( A-l) 


so that 


and 


r x r = 0 


r x r = L 


( A-2) 


where L is a vectorial integration constant. 
Consider, now, the result, 


_d / r_ 
dt \ r 


r‘ r - r\ r • r 


( r • T ) __ (~r • r ) f - r (r • r) r x(r • r) 


( A — 3) 


Combining equations (A-2) and (A -3) gives 


_d_ 

dt 



- r x L 


( A-4) 


54 



Similarly, we compute 


d /r*- -7* N -r: t- rxL , . 

— (rx L ) = r x L = - /x — ^5— • (A- 5 ) 

Multiplying equation (A- 4 ) by p and subtracting the result from (A- 5 ) gives 
d / r \ d ,-r -► . 

" STVT)" TTt (r x L> 

yielding 

r* x L = w ( A-6) 

where W is a second vectorial integration constant. 

It is instructive to identify W at this time. From (A-6) 

r • (rxL) = pr+r*W- 

But 


r • (rxL) = L(rxr) = L* L = L 2 

so that 

ju r + T • W = L 2 . 

Let ip be the angle between the vectors r and W so that 


L 2 = r ( n + w cos ip ) 
which gives 


JL 2 

r = * 

w 

1 + — cos ib 


(A- 7 ) 


(A-8) 
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The expression for the radius in the two-body problem described by 
spherical coordinate angles (p and 0 is [ 2] 


a(l - e 2 ) 

1 + e cos 0 cos 4> 


( A-9) 


Comparing equations (A-8) and (A-9) shows that 


w = jue (A-10) 

where e is the eccentricity of the orbit. In the standard choice of coordinates 
where perigee corresponds to ip = 0 we find that W is a vector of magnitude 
lie pointing toward perigee. 

Equation (A-6) now yields the z component of the vector W which, in 
parabolic coordinates, becomes 


(W »z = 


U a 


m 2 (|+rj) 


( p , 2 - »/)■ 


(£ -n) p 




2m ^ 


L 2 


^ ~ 1 
£+7] 


(A-U) 


Ihe right side of equation (A- 11) can be seen to differ from the portion 
of equation (51) which is in curly brackets only by a factor of m 2 . (Notice that 
1 = L if A = 0. ) The identification of /? is then complete. It is a modified form 
of one component of Laplace's integration constant just as E and 1 were, 
respectively, modified forms of the energy and angular momentum. 
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